home *** CD-ROM | disk | FTP | other *** search
/ Libris Britannia 4 / science library(b).zip / science library(b) / MATH / LGSOLV1.ZIP / README.DOC < prev   
Text File  |  1993-10-03  |  8KB  |  210 lines
  1.         LGSOLV - a set of FORTRAN-functions for
  2.              Large Linear System solving.
  3.         =======================================
  4. Authors: P.G.Akishin, A.P.Sapoznikov.  LCTA,JINR, Dubna,Russia.
  5. E-mail:     akishin@main1.jinr.dubna.su
  6.          sapoznikov@main1.jinr.dubna.su
  7.  
  8.  
  9.         If you have to solve Large Linear Systems, then a
  10. program described below will be useful for you. If you
  11. haven't - you simply will get information about a typical Tool
  12. for such problems.
  13.         The main problem for programmers working witn large
  14. matrices is how to keep them in computer's memory. Being
  15. written in a disk-storage, matrix involves a lot of disk
  16. exchanges, because all Matrix-Algorithms use it "row-by-row"
  17. and "column-by-column" at the same time.  Here we're proposing
  18. an original modification of known Gauss method for pivot element
  19. elimination that uses internal rows-transposition-vector and
  20. works with matrix strictly "col-by-col".
  21.         Our algorithm factorizes the source matrix so that
  22. created factor-matrix may be used several times for quick
  23. solution of Linear System with any Right Parts.
  24.  
  25.       Your Main Program ought to do only two things :
  26.   1.  to prepare functions calculating matrix elements
  27.       and system right parts;
  28.   2.  to reserve 2 buffers in memory, the larger they are,
  29.       the less will be a number of exchanges.
  30.  
  31. So, let us start.  The system  A*X = C  is being solved, where
  32. A is a square matrix [ NDIM : NDIM ],   X,C - are
  33. NDIM-vectors.  In practice it's often necessary to solve the
  34. whole family of NC systems:
  35.  
  36.                A * X(k) = C(k)          k=1,2,...,NC
  37.  
  38. with the same matrix A.  Historycally FORTRAN was being used
  39. for similar problems, that's way our algorithm is, of course,
  40. FORTRAN-written.  Here's our Operation Set for Large Linear
  41. System solving and results getting:
  42.  
  43.         INTEGER FUNCTION LGMINP( lbuf, mpf, ndim, elmatr )
  44.         INTEGER FUNCTION LGMINV( lbuf, mpf, ndim, elmatr )
  45.         INTEGER FUNCTION LGSOLV( nc, elvect )
  46.         REAL*8  FUNCTION G1SOL ( row, nrp )
  47.         REAL*8  FUNCTION G1RP  ( row, nrp )
  48.         REAL*8  FUNCTION G1EMAT( row, col )
  49.         REAL*8  FUNCTION G1FMAT( row, col )
  50.         INTEGER lbuf,mpf,ndim,nc,row,col,nrp
  51.         REAL*8  elmatr,elvect
  52.  
  53.    ndim - System Dimension
  54.    nc   - number of System Right Parts
  55.    lbuf - size of 2 working buffers declared in main program as
  56.              common/buf1/b1(lbuf)
  57.              common/buf2/b2(lbuf)
  58.              real*8 b1,b2
  59.           You are giving to SOLVER as many space, as you can,
  60.           but no less than  ndim.
  61.    mpf  - matrix packing factor (4 or 8).  All the calculations
  62.           use Double Precision but only mpf high bytes of each
  63.           word will be written on disk. Default=8. Usage mpf=4
  64.           requires less working disk-storage but decreases the
  65.           result's accuracy.
  66.    elmatr - user function giving a Source Matrix element
  67.    elvect - user function giving an element of Right Part
  68.  
  69.    LGMINP inputs a Source Matrix  A  invoking  elmatr(i,j)
  70.           ndim**2  times in a following sequence:
  71.           (1,1),(2,1),...,(ndim,1),(1,2),(2,2),...,(ndim,2)...
  72.           i.e. "col-by-col",  and writes it to Scratch-file
  73.           using logical number 1.  Than LGMINP prepares two
  74.           triangle matrices, L and R,  getting  A = ~L * R, and
  75.           writes them to Scratch-file using logical number 2.
  76.           Both files are created in the current user's working
  77.           directory automatically and are deleted when user
  78.           program ends. Buffers b1,b2 are used for the
  79.           disk-exchanges.
  80.    LGSOLV inputs all the System Right Parts  C  invoking
  81.           elvect(i,k)  ndim*nc times in a following sequence:
  82.           (1,1),(2,1),...,(ndim,1), ..., (ndim-1,nc),(ndim,nc)
  83.           and adds them to File-2.  Than LGSOLV computes the
  84.           System Solutions  X = ~R * ( L * C )  and adds them
  85.           to File-1.
  86.    G1SOL - these 4 functions are simply extractors of information
  87.    G1RP    that SOLVER's files hold.  G1SOL gives row-th element
  88.    G1EMAT  of nrp-th solution.  G1RP gives row-th element of nrp-th
  89.    G1FMAT  Right Part.  row,col=1,2,...,ndim;  nrp=1,2,...,nc;
  90.           G1EMAT gives (row,col)-th element of Source Matrix.
  91.           G1FMAT - the same for factor-matrix (we don't imagine
  92.           whom it will be useful,  it was made "for symmetry").
  93.    LGMINV realises a special problem : Matrix Inversion.
  94.           Having inputing the matrix, it just starts to solve
  95.           system  A * X = E  and writes all  ndim  solutions
  96.           to File-1 instead of Source Matrix.  In this case
  97.           G1EMAT gives elements of Inverse Matrix but
  98.           G1SOL,G1RP - are nonsensible.
  99.  
  100.  
  101.  
  102.  WARNING: as we hold information "col-by-col", be careful during
  103.     extraction : if you extract "row-by-row" then it takes a
  104.     disk-exchange for almost every extracted element !
  105.  
  106.  
  107.                     Return Codes
  108.  
  109.  All the operations don't print any message or stopped.
  110.  REAL*8 entries simply return  0.0  if their parameters are wrong.
  111.  INTEGER entries return codes are :
  112.         0 - all right
  113.         1 - Source Matrix is singular
  114.         2 - user buffer's size    lbuf < ndim
  115.  
  116.  
  117.                     An Example
  118.  
  119.    In your application program there ought be the following :
  120. C
  121. C  .....
  122. C
  123.       PARAMETER(lbuf=...)
  124.       COMMON /buf1/ b1(lbuf)
  125.       COMMON /buf2/ b2(lbuf)
  126.       EXTERNAL fun1,fun2
  127.       REAL*8 b1,b2,fun1,fun2,G1SOL,G1EMAT,G1RP
  128. C  .....
  129. C       1.  Input a System Dimension and Matrix into LGSOLV.
  130. C           here  ndim = <your System Dimension>,
  131. C                  mpf = <Matrix Packing Factor (4 or 8) >
  132. C
  133.       ner = LGMINP( lbuf, mpf, ndim, fun1 )
  134. C
  135. C  .....
  136. C       2.  System Solving for NC given Right Parts :
  137. C
  138.       ner = LGSOLV( nc, fun2 )
  139. C
  140. C  .....
  141. C       3.  Results getting : i,j=1,2,...,ndim;  k=1,2,...,nc
  142. C           (i-th element of k-th solution, k-th Right Part
  143. C            or j-th  Source Matrix column ) :
  144.       x = G1SOL( i, k )
  145.       c = G1RP ( i, k )
  146.       a = G1EMAT( i, j )
  147. C  .....
  148.       END
  149.  
  150.       REAL*8 FUNCTION fun1( i, j )
  151.       INTEGER*2 i,j
  152. C       This function gives A(i,j) - element of System Matrix
  153.       fun1 = ...
  154.       RETURN
  155.       END
  156.  
  157.       REAL*8 FUNCTION fun2( i, k )
  158.       INTEGER*2 i,k
  159. C       This function gives C(i,k) - i-th element of
  160. C       k-th  System Right Part Vector
  161.       fun2 = ...
  162.       RETURN
  163.       END
  164.  
  165.  
  166.                     The Demo Test
  167.  
  168.       Our Demo Test is prepared for IBM PC.
  169.       When you run  DEMOTEST.EXE  file, you'll be asked about
  170.       current lbuf,mpf,ndim,nc  values and matrix type.
  171.       Test solves the system  A*X=C  choosing C corresponing
  172.       to known solution X(i)=i.
  173.       During the solving process you can see the dynamic of
  174.       source and factor matrices usage.  After process' finish
  175.       elapsed time, exchanges's counter and errors will be shown.
  176.       Errors are the difference between computed and known
  177.       solutions.
  178.       Using the 286-processor we didn't want to compile in LARGE
  179.       mode. Because of that you can't use  lbuf > 8000
  180.       In the other hand, it seems that 8000 will be quite
  181.       enough for all cases.  You can to experiment with lbuf
  182.       working with Demo Test. If you are going to do it,
  183.       don't forget to experiment with Matrix Packing Factor
  184.       also. You'll have seen, how the errors depend on MPF.
  185.       You imagine, of course, that the sensible NDIM for test
  186.       is about 50-100.  NDIM=400 requires about 250 sec. for
  187.       386 processor.
  188.  
  189.  
  190.                         LGSOLV Packing List
  191.  
  192.       1. DEMOTEST.EXE  - run file for IBM PC
  193.       2. DEMOTEST.FOR  - Source FORTRAN text of DEMOTEST Program
  194.       3. LGSOLV.OBJ    - compiled by MS-FORTRAN 5.00   for 286 processor
  195.       4. README.DOC    - this document
  196.  
  197.   If you want to join LGSOLV.OBJ to your own main program,
  198.   add the following empty subroutine:
  199.  
  200.       subroutine vis0(m1,m2,job)
  201.       integer*2 m1,m2
  202.       character*(*) title
  203.       entry vist(title)
  204.       entry vis(m1,m2,job)
  205.       return
  206.       end
  207.  
  208.   because LGSOLV Demo Version invokes these 3 subroutines
  209.   for solving process' visualisation.
  210.